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Abstract 

Analytical (rational) mechanics is the mathematical structure of Newtonian 
deterministic dynamics developed by D'Alembert, Langrange, Hamilton, Jacobi, 
and many other luminaries of applied mathematics. Diffusion as a stochastic pro- 
cess of an overdamped individual particle immersed in a fluid, initiated by Ein- 
stein, Smoluchowski, Langevin and Wiener, has no momentum since its path is 
nowhere differentiable. In this exposition, we illustrate how analytical mechan- 
ics arises in stochastic dynamics from a randomly perturbed ordinary differential 
equation dXt = b(X t )dt + edWt where Wt is a Brownian motion. In the limit 
of vanishingly small e, the solution to the stochastic differential equation other 
than x — b(x) are all rare events. However, conditioned on an occurence of 
such an event, the most probable trajectory of the stochastic motion is the solu- 
tion to Lagrangian mechanics with C — \\q — b(g)|| 2 /4 and Hamiltonian equa- 
tions with H(p,q) = ||p|| 2 + b(q) ■ p. Hamiltonian conservation law implies 
that the most probable trajectory for a "rare" event has a uniform "excess ki- 
netic energy" along its path. Rare events can also be characterized by the prin- 
ciple of large deviations which expresses the probability density function for Xt 
as f(x,t) = e~ u ^ x,t ^ e , where u{x,t) is called a large-deviation rate function 
which satisfies the corresponding Hamilton- Jacobi equation. An irreversible dif- 
fusion process with Vxt/O corresponds to a Newtonian system with a Lorentz 
force q = (V X 6) X q + |V||6|| 2 . The connection between stochastic motion 
and analytical mechanics can be explored in terms of various techniques of ap- 
plied mathematics, for example, singular perturbations, viscosity solutions, and 
integrable systems. 

1 Introduction 

Dynamics as an analytical concept is one of the most important contributions of math- 
ematics to modern thinking. Currently there are three fundamentally different types of 
dynamics: classical, quantum, and stochastic [1|. Quantum dynamics aside, classical 
dynamics in term of deterministic "trajectories" of a system, continuous or not, is still 
the dominant model in quantitative science and engineering. However, rapid develop- 
ment of quantitative biology in recent years, and the sustained interests in statistical 
physics, has pushed the third, stochastic dynamics, to the forefront of applied mathe- 
matics. For recent reviews on Darwinian dynamics and the Delbriick-Gillespie process 
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for cellular biochemical systems, see (HE). Novel and even non-orthodox stochastic 
dynamic approaches to quantum phenomena can be found in [3]|4][5). 

Describing the celestial mechanical system of a few interacting bodies, Newton's 
equations of motion in terms of classical dynamics is one of the most successful math- 
ematical models known to mankind. An in-depth study of the subject in terms of an- 
alytical mechanics exposes one to a wide range of applied mathematical theories and 
techniques: Hamiltonian systems, Lagrangian principle, and Hamilton-Jacobi equa- 
tions are several shining jewels of the treasure box ||6]|2l[8]. 

In this exposition, we shall illustrate that much of these classical, deterministic 
mathematics also emerge naturally in the theory of stochastic dynamics, when one is 
interested in the relationship between stochastic and deterministic dynamics. This is 
reminiscent of the semi-classical theory of quantum dynamics developed in the 1970s 
|9|. More importantly, not only are they useful as tools for solving problems, several 
quantities also acquired a strong probabilistic meaning, such as large deviation rate 
functions and most probable paths. 

Even though stochastic dynamics is often described in terms of its probability dis- 
tribution changing with time, we need to emphasize, at the onset, that neither the distri- 
bution perspective, nor a trajectory perspective, is a complete description of a stochastic 
process. The mathematical notion of a random variable, developed by A. Kolmogorov, 
can not rest solely on its distribution function, nor its realizations 11101 . It is truly an 
independent new mathematical object with deep philosophical consequences^] 

2 Stochastic Dynamics in Terms of Stochastic Differen- 
tial Equations and Diffusion Processes 

One of the extensively studied problems that connect stochastic and deterministic dy- 
namics is an ordinary differential equation with small random perturbations ifTTI . Let 
us consider a diffusion process X t with the stochastic differential equation 

dX t = b(X t )dt + V2~edB t , (1) 

and the corresponding Kolmogorov forward equation 

In the limit of e = 0, |2]) is formally reduced to the first-order partial differential 
equation (PDE): 

«fM_-£ W .)/(., t >), (3) 



1 Two issues immediately come to mind: (i) The mathematical theory of probability requires all possible 
outcomes being known a priori, in the very definition of a random variable. This makes the concept of a 
random variable only retrospective. In statistics, this is intimately related to the concept of a prior distri- 
bution; and in economics this distinguishes risk from uncertainty, (ii) Classical dynamics has trajectories 
but only singular distributions; quantum dynamics has distributions but no trajectories due to Heisenberg's 
uncertainty principle; stochastic dynamics requires both perspectives. 
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which is equivalent to, according to Liouville's theorem, a nonlinear ordinary differen 
tial equation (ODE) 

dx(t) 



dt 



= b(x(t)). 



(4) 



It is known that on any finite time interval t g [0, T], the solution to Eq. ([T} with initial 
value Xq — x Q approaches the solution of the ODE with probability 1 Bill . 
If we write f(x, t) = e™^*', then the linear PDE ^ becomes 



dw(x,t) dw(x,t) db(x) 



01 



dx 



dx 



(5) 



a linear, first-orderpartial differential equation. 

Now for Eq. let us assume f(x, t) — e -««(*.*)/' 
ansatz. Then we have 

du e (x, t) 



dt 



du e 
dx 



b(x) 



du e 
dx 



, widely known as the WKB 

d 2 u e db(x 



dx 2 



dx 



(6) 



Therefore, to leading order, if lim e 
one has 

du(x, t) 
dt = 



du(x, t) 
dx 



u(x, t) exists and is differentiable, 

du(x, t) 



-b{x) 



dx 



(7) 



Eq. (j7]i is widely called the Hamilton- Jacobi, or Eikonal, equation (HJE) lfT2l[T3l . Note 
that this equation is different from Eq. (|5j: It is a nonlinear, first-order PDE. As we 
shall discuss below, the solution to |5]) is the limit of equation |2]) when e — >• 0. The 
solution to (|7| is the convergence rate of that limiting process. u(x,t) is called the 
large-deviation rate function. 

To relate a continuous random variable X with probability density function fx(x) 
to a deterministic quantity, the expected value E[X] and modal value x* , with /(x*) > 
f(x), are often taken as the counterpart. Let us now consider a local minimum of the 
function u(x,t) located at x = x*(t). How does the location and the value of the 
minimum change with time? According to Eq. 



du{x* (t),t) 
dt 



du{x, t) 
~lft 

/ du(x, t) 
~ V dx 



du(x, t) 
dx 



b(x) 



dx*{t) 
dt 



x=x* (t) 



du(x,t) ( du(x,t)\ dx*(t) 



dx 



dx 



dt 



(8) 



X — X* (t) 



= 0. 

Also from Eq. |7]) we have: 

d fdu{x*{t),t) 



= 



dt 



dx 



d 2 u{x,t) d 2 u(x,t) ( dx*{t) 



dxdt 



dx 2 



dt 



(9) 
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Therefore, 



dx*(t) 
dt 



d 2 u(x,t)\ d 2 u(x,i) 



dx 2 



d 2 u(x, t) 
dx 2 

d 2 u(x,t) 



dxdt 

- x—x*{t) 

du(x, t) d 2 u(x, t) 



b(x) 



b(x*(t)). 



dx 2 



dx 



x—x* (t) 



dx 2 



db(x) du(x, t) 
dx dx 



(10) 



So indeed, the modal values follow the ODE Q. For an ODE with multiple domains 
of attraction, they correspond to a multi-modal distribution. Furthermore, 



d_ ( d 2 u(x*(t),t) 
dt [ dx 2 



v d 3 u(x,t) d 3 u(x,t) fdx*(t) 



dx 2 dt 



du 

dx 



dx 3 \ dt 
d 3 u d 2 b(x) 



r*(t) 



dx 3 



dx 2 



db{x)d 2 u (dx*{t) 



dx dx 2 



dt 



b(x. 



d 2 u 
dx 2 

d 3 u 
dx 2 



= -2 



fdb(x*(t)) d 2 u{x*(t),t)\ d 2 u{x*(t),t) 



\ dx 



dx 2 



dx 2 



.(11) 



Hence, the value u(x* (t), t) does not change, and its location follows the ordinary dif- 
ferential equation x — b(x). A local minimum of the function u(x, t) follows the cor- 
responding deterministic ODE. Furthermore, the curvature about it follows Eq. 
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Since u(x, t) is the rate of convergence of a normalized probability distribution for 
X\ when e — > 0, it is a non-negative function with its minima necessarily zero. In 
fact, except for the very critical condition known as phase transition, the minimum is 
unique. Eqs. <JsJ and ( 10 1 state that if an initial u(x, 0) is non-negative with a global 
minimum zero, u(x, t) will remain non-negative with global minimum zero. In other 
words, the properties of being a large-deviation rate function are preserved. 

Finally, observing that the values of local minima are related to the probability 
associated with each "attractor", Eq. ([8]) states that diffusion processes in different 
"attractors" are almost reducible in the limit of e — > 0. Generically speaking, besides 
the dominant attractor with the global minimum of u(x, t), the probability of each 
attractor with a local minimum vanishes. This results in the Law of Large Numbers. 
However, conditioned upon being outside the dominant attractor, there will be another 
global attractor. These states are known as metastable in statistical physics. 
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3 The Probabilistic Interpretation of u(x, t) 

We now give the precise meaning for u(x, t). Let us denote the solution to Eq. (|2j, with 
initial condition f(x, 0) = 5(x — x a ). As the solution to the transition probability 
/ c (a;, <:jrc ) is the fundamental solution to the linear PDE. For fixed x, t, and x a , when 
e — > 0, one has 

lim/ e (a:,f) = <J (ic-xt), (12) 

e-s-0 

where a; t is the solution to the ODE x t = b(x t ) with initial value x Q . Or in a more 
authentic probabilistic notation: 



limPr{|X t (e;z )-^l > e} = 0. 



We therefore can introduce the rate of convergence: 

Pr{x < X t (e;x ) < x + dx} 



lim — e In 

e-S-0 



dx 



= u{x,t;x ). 



(13) 



(14) 



3.1 Laws of large numbers, central limit theorem and theory of 
large deviations 

Let us consider a sequence of iid (independent, identically distributed) Xj and their 
mean value 



1 

Z n = — / Xi 

n — ' 



(15) 



The weak law of large numbers (LLN) states that 

lim Z n =fj,*, n* = E[X\. (16) 

n—^oo 

That is, the probability density function for a continuous random variable 

lim f Zn {z) =8{z-n*), n* = E[X\. (17) 

n— > oo 

Furthermore, the central limit theorem states that 



lim —j=Sz n 



1 



Var[X}. (18) 



We note that the result in Eq. ( 18 i implies Eq. ( 17 1, which implies Eq. ( 16 1. 
The large-deviation rate function for the LLN given in Eq. ( 16 1 is defined as 

u(x) = lim — In fz.(x). 



(19) 



Note that the convergence of the probability distribution functions in Eq. ( 17 1 is always 
non-uniform. Asymptotics beyond all orders, therefore, necessarily enters the theory 
of large deviations. 
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According to Chernoff's formula [14), u(x) is related to the cumulant generating 
function (CGF), X(9), for the random variable X: 

\(8) =\nE[e- dx ] , (20) 

via the Fenchel-Legendre transform 

A(0) = sup{6»a; - u(x)}. (21) 

X 

Furthermore, we know that the Fenchel-Legendre transform of X(0), 

u*{x) = sup{x9 - A(0)}, (22) 
e 

is the affine regularization of u(x). u*(x) is a convex function; it is identical to 
u(x) in the neighbourhood of the global minimum of u(x). Cumulant (or the Thiele 
semi-invariants) expansion has recently found applications in renormalization-group 
approach to singular perturbation [ 15 1. 



From Eq. ( 22 1 and the basic properties of the Fenchel-Legendre transform, we have 

mmu*{x) = -A(0) = 0, (23) 

X 

„. (a)=o = A'(0)| e=o = /A (24) 

_ , = (^'WLor 1 = (VarM)" 1 = \. (25) 



x 



X — fjb 

Therefore, we have in the neighbourhood of global minimum of u(x), x = fj,*; 

u{x) = u*(x) = {X ~£ )2 + a 3 (x - (i*) 3 + -... (26) 

While u*{x) is convex, u(x) need not be. Hence, away from the global minimum 
x = /!*, u(x) can have many local minima. Let us denote their locations as xi and 
with corresponding local expansions 

u(x) =u{x t )+ a ^-{x-x t f + ■■■ , e=l,2,--- (27) 

in which u(x^j > and > 0. The fz n (x) then has an asymptotic expansion 

-nu(x)+o(n) / _ —nu(x)-\-o{n) 



f z f x ) = e""W w = J e -n«W+»W, (28) 

2na 2 
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Eq. (28 i implies Eq. ( 18 1 in the following sense: 



fzjx) 



e 2^ 



27TCT 2 



dx 



2ira 2 



-nu{x)+o{n) 



-e 2<r 2 



72 n (x- M ')^ 

-e 2,2 



27TCT 2 



27TCT 2 

-na 3 ((a;— p*) 3 +o(ri) _ 



-+o 



•vAi \ \/n 



1 







1 



That is, the convergence of the distribution in Eq. ( pL8) > is in Li and it is of the order 



Note that the non-convex parts of u(x), as given in Eq. (27i, only contribute to 
terms on the order of 



e -nu(x e ) / e -^JL( x - xe ) dx = 



(29) 



It is exponentially small, i.e., beyond all orders. 



3.2 Linear dynamics, Gaussian processes, and an exactly solvable 
HJE 



We now consider the SDE in ([T]) with linear drift, i.e., a Gaussian process: 

dX t = -bX t dt + V2edB t . 
The corresponding Kolmogorov forward equation is 

df(x,t) d 2 f , d 



dt 



and the WKB ansatz leads to the PDE 



du e (x, t) 
dt 



du e 
dx 



bx 



du e 
dx 





~d 2 u e ' 


+ e 


9a; 2 



(30) 



(31) 



(32) 
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for u e (x, t). We can expresses 



u e (x,t) = a(t) + {X 2(J ^\ 03) 



then we have a set of nonlinear ODEs: 

d 



a 2 (t) = 2(l-6er 2 (i)) , (34a) 



dt 
dn{t) 
dt 



-bfi(t), (34b) 



dt \<J 2 



Their explicit solution is 



a\t) = \+(° 2 ^)-^)^ 2b \ (35a) 
H{t) = n{Q)e- h \ (35b) 

«(*) = | ln (^|y) + a (°)- ( 35c ) 

Assembling these together, we have 

e-<**»' = . ^ exp{- ( V^; )2 | , (36) 

where A = ^/ 2ne<7 2 (0) e~ a (°'/ c is a constant. 

We note that in Eq. ( |34) , the e term only contributes to a(t). This is a very sur- 
prising result: The diffusive behavior characterized by a 1 (t) actually is determined by 
a first-order nonlinear PDE (|7|! As we shall see below, the solution to Eq. |7]l can be 
constructed via the method of characteristics. 



3.3 Large-deviation rate function on a circle 

We now consider dynamics on a circle 6 £ S 1 [0, 1] with periodic angular velocity 
6 = b(0) = b(8 + 1). In general, a saddle-node bifurcation on the cycle gives rise 
to a counter-clockwise or clockwise cyclic motion when either b(6) > or b(6) < 
W € [0, 1] |16|. We have shown above that the corresponding large-deviation rate 
function u(6) has its minima and maxima corresponding to the stable and unstable 
fixed points. We now illustrate that corresponding to this cyclic motion, u(9) becomes 
a constant on the entire circle. 

Without losing generality, we shall assume there is only one minimum and one 
maximum of u(x, t) and let x(t) and y(t) be the location of the minimum of u(9, t) 



9 



and the corresponding curvature. We are interested in how x and y behave in the 
infinitely long time limit. Following Eqs. ( 10 1 and ( 1 1 1 we have 

(IT 

— = b(x), (37a) 
d j- t = -2{b\x)+v)y. (37b) 

Fig. [TJgraphically shows the occurence of a Hopf bifurcation in this autonomous planar 
system in R 2 , corresponding to a saddle-node bifurcation of b{&) on the circle § . 

We now show for the case of x — b(x) > periodic solution as shown in Fig. [Tj3, 
with y(t) — > asymptotically. Since x(t) is a periodic function of t, b'(x(t)) is also 
periodic. Furthermore, 

Mt) = f b' (x(s)) ds = [ b '^ S }} ( ^jfT] ds = In b(x(t)) + const. (38) 
Jo Jo °K S )) \ ds J 

is also a periodic function of time. Then 



This equation can be re-written as 

d 
dt 



(y(t)e 2 ^) = -2 (y(i)e 2 *)) 2 e- 2 *>. (40) 



The right-hand-side of Eq. (gOjl is < and it is iff y{t) = 0. Therefore, y(t) -> if 
2/(0) > 0. In other words, the curvature approaches zero. The large deviation function 
along a limit cycle is a constant ifTTl . 



4 A "Fictitious" Classical Newtonian Motion 

Eq. (|7|l has the form of a Hamilton-Jacobi equation (HJE) for a classical Newtonian 
motion according to analytical mechanics OlIllHl- 

4.1 A Hamiltonian system 

In analytical mechanics, the solution u(x, t) to the HJE, a nonlinear PDE like Q, is 
called a principal function. u(x, t) furnishes the entire family of orbits, i.e., a flow, 
corresponding to an associated Hamiltonian dynamical system in terms of a system of 
nonlinear ODEs. 

The Hamiltonian system, with time-independent Hamiltonian^] 

H(q,p)=p 2 + b(q)p, (41) 

2 The Hamiltonian associated with a discrete birth-and-death process (B&DP) with birth and death rates 
u n anAw n is H(p,q) = u(q)e p + w(q)e~ p — u(q) — w(q) where u(q) andto(g) are the continuous limits 
of u„ and w„. Note this H (p, q) is reduced to Eq. {41) if p <g 1. The B&DP and diffusion corresponding 
to very different forms of the HJE leads to the diffusion's dilemma 1181 . 
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is a system of autonomous ODEs 



% = 2p + b{q), 



dH 



db(g) 



(42) 



It is a re-formulation of a Newton's "equation of motion" 



d 2 q ,db{q) 



dq dq 



(43) 



This "fictitious" analytical mechanical system associated with stochastic dynamics 
has been identified and explored by Graham and Tel lfT9ll . and in the even earlier work 
of Martin-Siggia-Rose ll20l . 

According to analytical mechanics, the Hamiltonian system also has a correspond- 
ing "Lagrangian" 



with the action functional, 



S [q(t);(t ,q(t ))^(t,q(t))} 



dq[s) 
dt 



b(q(s)) 



ds . 



(44) 



(45) 



In classical mechanics, the action has never had a "meaning" more than being a math- 
ematical device which provides trajectories for a mechanical motion. In the present 
work, however, e~ s °' e turns out to be exactly the probability of a path of the stochas- 
tic dynamics in ([T]i with vanishing e. For finite e, the generalized Onsager-Machlup 
functional is S e = S„ + (e/2)b'(q) ED [22] [23] |24] . The action functional in Eq. ([45 1 



plays a central role in Freidlin-Wentzell's theory of large deviations ifTTl . Also see 
||251 for a rigorous mathematical treatment and ll26ll for a monograph with an applied 
mathematical bend. 

One should not confuse the present Hamiltonian system with a more involved 
Hamiltonian formulation associated with stochastic dynamics recently proposed in 

ma. 



4.2 Path integral formalism for the probability of a diffusion pro- 
cess 



We now show heuristically that the path integral in Eq. (45 i represents the probability 
of a diffusion process according to Eq. Q. Statistical physicists have long used it 
as a useful mathematical tool, starting with Onsager and Machlup 12T1 for Gaussian 
processes and later by Haken, Hunt and Ross [22 , 24 1 for the general diffusion process. 
See EH1 for its application in the protein folding problem and [29 1 for a very recent 
study on nonequilibrium steady-state physics using this formalism. 
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First, let us consider the solution to a simple one-dimensional diffusion process 
defined by Eq. ([I]) with a linear drift bx. The corresponding Fokker-Planck equation 
|2} has the fundamental solution 



f(x,t+ At\x',t) = 



AneAt 



( e 2b&t -l \ 
^ 2fcAt J 



■ exp 



(x-x'e bAt ) 2 
~2(e/b) (e 2bAt - 1) 



: exp 



(or - x' - bx'At) 



2 1 



AeAt 



(46) 



^ Aire At (1 + bAt) 

Now for the solution to the general Eq. (|2| with drift b(x), the probability density of a 
trajectory is approximately 

f (x n ,nAt;x n -x, (n - l)At; • • • ,xi,At\x o ,0) 



1 

\/47reAi 



n 



exp 



(xi - Xi-i - 6(a;j_i)At) 2 
4eAi 



•ln(l + 6'(xi)At) 



.4 exp < — — , 
I 4e J Q 



(dx{t) 
\ dt 



b{x{t))\ +2eb'(x(t)) 



dt 



(47) 

where t = nAt and A is an appropriate normalization factor. Note that according to 
Wiener's theory of diffusion, the stochastic trajectory x(t) is nowhere differentiable 
QUI . Hence, the use of dx(t)/dt in Eq. |47| i is heuristic, and only becomes mathemat- 
ically meaningful in the limit e = 0. See 1 25 1 and [11] for a rigorous treatment. 



4.3 Conditional probability interpretation of "excess kinetic en- 
ergy" 

A Hamiltonian system has a conserved quantity: the H(p(t), q(t)), along each and ev- 



ery trajectory. If we re-arrange the first equation in (42 1, we obtain p = (q — b(q)) /2. 
Then, 

H{p, q)=p(p + b(q)) = \<A- b{q)) (q + b(q)) = \ {q 2 - b 2 (q)} , 



that is 



q 2 = b 2 (q)+4H, 



(48) 



in which H is constant along a trajactory. Comparing this with the noiseless trajectory 
q = b(q) with H = 0, we see there is a constant 4H added to the q 2 . This result can 
be interpreted as follows. 

We shall call the square of the velocity q 2 the "kinetic energy". A noiseless trajec- 
tory, i.e., the one with H = 0, follows the differential equation q = b(q). This means 
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that with starting time t = t\ at position qy, the noiseless trajectory will be precisely at 
q 2 at time t%. Any other trajectories arriving at q 2 with a different time are impossible 
for deterministic dynamics, and are rare events when e is small. If, however, one ob- 
serves such an rare event: motion from q\ to q 2 with a time t 2 ^ t 2 , what will be the 
"most probable trajectory among all possible q(t) with q{t\) = q\ and q (t 2 ) = q 2 n - 
This is a problem of conditional probability. Among all the rare trajectories, 

{<z(*)l(K*i) = 1lA (*2) = <llM # *2>, 

the most probable one follows a solution to the Hamiltonian system with an appropriate 
H^0\ 

More explicitly, since q(t) arrives at q 2 with a time different from t 2 , say t 2 < t 2 , 
then it has to be speeded up compared to q(t). On the other hand, if t 2 > t 2 , it has to be 
slowed down. The solution to the Hamiltonian equation states that the most probable 
trajectory is the one having a constant amount of excess kinetic energy (AH) added or 
reduced along the trajectory. In the theory of probability, this is a consequence of the 
van Campenhout-Cover theorem PP . 

The fictitious Hamiltonian system, therefore, is the consequence of stochastic dy- 
namics conditioned on the occurence of a rare event. The kinetic energy exists "retro- 
spectively" in some stochastic dynamics. 



5 Solutions to the Hamilton- Jacobi Equation (HJE) 

Since the HJE is a nonlinear PDE, there is no systematic method to obtain its solution 
for arbitary initial data. Rather, there are classes of solutions one can obtain; exactly or 
approximately. Since this is a rather developed area of applied mathematics, we shall 
only touch upon some issues highly relevant to our mission. 



5.1 A class of exact solutions 

The standard way to solve a HJE is precisely by solving its corresponding Hamiltonian 
dynamics in terms of the ODEs. Since the Hamiltonian H(q,p) does not depend ex- 
plicitly on time, one can verify that u(x, t) = uq(x) — Et is a solution to HJE (|7| if 
Mo (x) is a solution to 

H (duoix) jX \ E _ 



d:r 



Therefore, 



u(x, t)= [ p(q)dq - Et = - f (-b(z) ± y / b 2 (z)+AE\ dz - Et. (50) 



Initial data uq(x) satisfying Eq. (49 1 is called characteristic initial data 

For E = 0, u(x, t) = 1*0 (x) Vt. Hence it is a stationary solution to the HJE. In this 
case, Eq. d49|) yields 

r 

u st {x) = - / b(z)dz (51) 
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and a trivial solution uq(x). The solution in Eq (51 1 corresponding to E = is the 
expected stationary solution. 



5.2 Solutions via characteristics 



The solutions given in Eq. (50 1 is in a special class. Certainly not any initial data 



u(x, 0) is characteristic. More importantly, the u(x, t) in (50 1 can not satisfy the basic 
properties of a large deviation rate function: It has to be non-negative with its minimum 
exactly being zero. Therefore, we need to look for other possible solutions to the HJE 
|7} corresponding to noncharacteristic initial data [32|. 
According to Evans ll32ll . let 

du du 

— i y y, — i — ^ u M> u, x i ^ x. 

at ox 

The the HJE |7) is in the form of 



F{y 1 z, u,t,x) = y + z 2 + zb(x) = 0, 



with 



y(s 

i(s 
x(s 



= -FL = -z 



db(x) 
dx 



= F' y y + F' z z = y + (2z + b(x)) z, 

= Fy(y,z,t,x) = 1, 

= F' z (y,z,t,x) = 2z + b(x). 



(52a) 
(52b) 
(52c) 
(52d) 
(52e) 



Note that y = —(z 2 + zb(x)) is the Hamiltonian —H(x, z). The pair of equations 



(52 1) and (52 ;) are the Hamiltonian system in Eq. (42 1. u(s) = xz — H(x, z) is the 



corresponding "Lagrangian". 



5.3 Phase portrait of characteristic lines of Hamiltonian system 

In our case, Eq. (42 1 can be solved since p 2 + b(q)p = E is a constant of motion: 



dq 
~dt 



±y/b 2 {q)+4E, 



i.e., 



dq 



Vb 2 (q)+4E 



±t + const. 



(53) 



When E = 0, this is precisely the solution to the ODE dq/dt — ±b(q). Fig. |5]shows 
the phase portrait of the Hamiltonian dynamics. 

The phase portrait in term of p as a function of q is 



p(q) 



-b(q) ± ^V(q)+±E 



(54) 
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5.4 Solution to HJE with noncharacteristic initial value 

As a large-deviation rate function, one appropriate initial condition for u(x, t) should 
be u(x, 0) = 0. It is clear that u(x, t) = is a solution to the HJE, but this is not a 
meaningful one. Therefore, one could be interested in the solution to the HJE with an 
infinitesimal initial data u(x, 0). Therefore, initially, one can linearize Eq. d6l: 



duo(x,t) du 



d 2 u db(x) 
dx 2 dx 



0(e 2 ). (55) 



The exact result for a Gaussian process can provide some insights: Note that if 
cr 2 (0) = oo, then the initial u(x, 0) = 0. 

5.5 HJE on a circle 

We again consider 9 € S 1 [0, 1] with periodic angular velocity 9 = b(9). The b{9) on a 
circle can be decomposed as 

with a differentiable periodic potential U{9) = U(9 + 1), 9 e S 1 , The corresponding 
Hamiltonian equation on a torus is 



= 2u + b(9) = 2io + b - d ^ : 
db{9) _ d 2 U{9) 



(57) 



with the Hamiltonian H(6,u) = uj 2 + ub -ljU'(0); {9,uj) e S 2 . 

In general, there are two types of fixed points (6, uj)* in system (57 1: to* = and 
U'{9*) = b , or U"{9*) = and to* = (1/2)(17'(0*) - 6 ). These two types are the 
same for systems with b — (i.e., a gradient system on S 1 ). For a rotational system 
with b(9) > or < on the entire S 1 , the first type does not exist. Therefore, taking 
the index theory for a C 1 -vector in a plane into consideration ll33l . the occurence of a 
cyclic motion in 9 = b{9) with changing bo corresponds to an annihilation of two fixed 
points of the two types. 

6 High Dimensional Cases: Momentum, Entropy Pro- 
duction and Nonequilibrium 

From a standpoint of the theory of Markov processes, the one-dimensional system in 
([T]i can only reach a time-reversible stationary process |34|. This result corresponds to 
the statement that an ODE on M. 1 is always a gradient system: x = b(x) = —dU (x)/dx 
with U(x) = — J b(x)dx. For autonomous ODE systems in higher dimensions, a limit 



cycle can occur. This corresponds to nonequilibrium phenomena in system (60i in 
dimension 2 or higher. 



15 



We now consider an TV-dimensional diffusion process with the stochastic differen- 
tial eqaution 

dX t = b(X t )dt + V2~edB tl X t £ R N , (58) 

and a coiTesponding Fokker-Planck equation for the probability density function f(x, t)dx 
= Pr{x < X t < x + dx}, 



d 



f(x,t) 



d 2 

e-^f(x,t) 



i) 



(b(x)f(x,t)). 



dt ' dx 2 ' dx 
6.1 HJE in two-dimensional systems 

For a two-dimensional SDE: 

dX t = b x (X t ,Y t )dt + V2edB^\ dY t = b y (X t .Y t )dt + V2edBl 2 \ 
the corresponding fictitious Hamiltonian is 

Py) = Pi +P 2 y + b x (qx,q y )Px + b y (q x ,q y )p y , 
so the Hamiltonian dynamical system is 



d 

dt q * = 


dH 

dp x 


= 2p x 


+ b x (q x ,q v ), 


d 

Jt q y = 


dH 

d Py 


= 2p y 


+ b y (q x ,q y ), 


d 


dH 

dq x 


= ~Px 


db x db y 
dq x Py dq x 


d 

dt Py = 


dH 

dq y 


= -Px 


db x db y 
dq y Py dq y 



Its corresponding Lagrangian is 

£ [qx,q y ,qx,q y } 
and the equations of motion are 

Cjx = 

q y = 



q x 




q y - by 



dby _ dbx 
dq x dq y 



db x 



dq x dq y 



In vector form Eq. ( 64 1 can be written as 

d 2 



^m = (v,b),^t) + 1 -vfbr. 



(59) 



(60) 



(61) 



(62) 



(63) 



(64) 



(65) 
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A Lorentz magnetic force like term arises if the vector field b(q x , q y ) is non-conservative. 

For a system of arbitary dimension, the vector form Eq. ( |65j ) is still valid if one 
interprets 

' dbi dbj\ dqj(t) 



dt. 



6.2 Force decomposition, momentum and entropy production 

We now show an interesting relation between the non-gradient b(x) and the large devi- 
ation rate function u st (x) for the stationary diffusion. 

We consider the case of only a single attractive domain with only one stable fixed 
point at 0. Suppose the TV-dimensional vector field b(x) admits an orthogonal decom- 
position |fl9l : 

b(x) = -VU(x) + £(x), (66) 

where the function U (x) is continuously differentiable, VU (x) ^ for x ^ 0, and 
the inner product in R , £(x) ■ VU(x) — 0. Then the large deviation rate function 
u st (x) = U(x)+const., and the unique extreme of the action functional So[x(t) : 
(to, 0) — > (T, x*)] is given by the equation 

^- = VU{x) + l{x), 

where t e [t Q ,T] and x(t ) = 0, x{T) = x* . 

Recall that u st (x) is the stationary solution to the Hamiltonian-Jacobi equation 

9U( ^ = -||V U (:M)|| 2 - Hx) ■ Vu(x,t), (67) 

associated with Hamiltonian H(q,p) = \\p\\ 2 + b(q) ■ p. 
Hence the Hamiltonian dynamics follow 



— = 2 Pi +bi(q) 
dpi 

dH db(q) 
~ is = -P- 



(68) 



dq t 

The action functional So is just the path integration of the corresponding Lagrangian 

C[ q ,q] = \p q -H( q ,p)]=J2 f4i ~ bi{qr ' 

i 

It is easy to calculate that along the classical trajectory associated with dx/dt = 
VU(x) + £(x), we have 

p(t) - -\7u st (x(t)); 

(69) 

q(t) = VU(q(t))+£(q(t)). 
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The Hamiltonian of this optimal dynamics is always zero, and the associated trajecto- 
ries cross all the fixed points of the deterministic dynamic system (ffl. 



In applied stochastic dynamics, the stationary solution u st (x) to Eq. ( 67 1 can be 
considered as a "landscape" for the dynamics |35l[T8l . Then Eq. ( |69"| l provides a very 
novel "meaning" for the conjugate momentum in the fictitious Hamiltonian system: 
It is the force associated with the landscape. Moreover, (p(t) , £(q(t)) = 0, i.e., the 
Lorentz force £(x) is perpendicular to the momentum. 

For a multi-dimensional diffusion process ( |58| l the entropy production rate for the 
stationary diffusion process is defined as 



3 p = - J \\b(x) — eV log7r(x; e)|j 2 7r(x; e)dx, 



(70) 



in which tt(x; e) is the probability density function for the stationary process. There- 
fore, e p = if and only if £(x) = 0. e p ^ is widely considered to be a fundamental 
property of a nonequilibrium steady state; e p = implies that the stationary process 
X t is time reversible |34|. 

When e tends to zero, if eVlog7r(x; e) converges, then the u st (x) is related to the 
stationary probability density through a Boltzmann-like relation 



limeVlog7r(a;;e) = -V/(i), 



(71) 



and then 



||6(x)-eVlog^(x;e)|| 2 ^ ||£(x)|| 2 , 



(72) 



the left-hand-side of which is inside the integral in Eq. ( |70| ). Therefore, asymptotically 
we have 

e p « 1 J \\£(x)\\ 2 e- u ^/ e dx. (73) 

Cases with multiple attractive domains are much more complex due to the problem 
of turning points and boundary layers in the limit of e tending to zero l36l l37l . 



6.3 Processes with time-reversal 



We now consider the diffusion process corresponding to Eqs. (58 i and (59 1 with a 
time-reversal. The corresponding forward equation for the reversed process is 



d ~ d 2 ~ d 



2e ^ ' - b(x) ) f(x,t) 



(74) 



where n(x) is the stationary density for Eq. (59 1. In other words, the corresponding 
drift is 

K*)=*^^-6(*). (75) 



In fact, both b(x) and b(x) can be written as 

91n7r(x) 



b(x) 



dx 



+ l(x), b(x) = e — ^-L-£( x ), 



(76) 
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with 

t{x) = Kx)-e™%&. (77) 
If we denote the stationary distribution 

tt(x) =e~ u{x)/ \ (78) 

then b(x) — ~VU(x) + £{x) and b(x) = —VU(x) — £(x). Their corresponding 
fictitious Hamiltonians are 

H(p,q)=p 2 +p(-VU(q)+£(x)) (79) 

and 

H(p,q)=p 2 +p(-VU(q)-e(q)). (80) 
The corresponding "equations of motion" are given by Eq. ( |65j ): 

^q(t) = ± (V X l{q)) x j t q{t) + \v (\\VU{q)\\ 2 + Mil 2 ) ■ (81) 
In the case of reversible processes, £(x) — and b{x) — b(x) — —VU(x). 

7 Interpretive Remarks 

Stochastic dynamics following trajectories defined by Eq. (|58|, with the time-evolution 



of the corresponding probability distribution characterized by Eq. ( 59 1, has a continu- 
ous but nowhere differentiable trajectory. Diffusion represents motions of a particle in 
a highly viscous medium in which inertia is completely lost instantaneously |[38l . At 
any given position x £ R^, the particle can move in any possible direction with any 
possible speed, but the mean velocity is (dx)/dt = b(x). There is no momentum in the 
classical sense. 

The analytical mechanical structure hidden in the stochastic dynamics discussed in 
the present review, however, is Newtonian. We note that the acceleration for determin- 
istic dynamics 

§-«.), (82) 



is 



d 2 Xi d ^ dbi(x) 

j=i j 

I 

1 



3=1 ' 3=1 



(V x b{x)) x b(x) + ^V||6(x)|| 2 . (83) 
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This is exactly Eq. ( 65 1. But how should one interpret the inertia and momentum in the 
fictitious Newtonian motion? How can we interepret the Hamiltonian as a conserved 
quantity in a stochastic trajectory? 

Mean dynamics and probability moment closure problem. The mean behavior of 
nonlinear stochastic dynamics like ( |58| , while being a deterministic function of time, 
can not be represented by a simple ordinary differential equation. For a stochastic 
process x(t): 

<L{x(t)) = {b(x))^b({x(t))). 

A traditional approach to resolve this problem is to introduce higher-order moments for 
the distribution of x(t), and to express the stochastic dynamics in terms of the mean, 
variance, third moments, etc. The method of moment closure was introduced to reduce 
this infinite hierarchy to approximately a finite system. 

Such an approach encounters significant difficulties if a nonlinear b(x) has multiple 
domains of attraction. In this case, the mean dynamics reflects two fundamentally 
different behaviors on very different time scales: the intra-attractor dynamics and inter- 
attractor dynamics |0. For a small e, the latter are rare events. 

Kinetic energy without momentum. The equations in ( 82 1 and ( 83 1 seem to suggest 
another line of hierachical characterization of stochastic dynamics. The Hamiltonian 
dynamics corresponding to Eq. ( 83 i yields Eq. < [82] > when H = 0. For all other dynam- 
ics with H ^ 0, they are impossible for Eq. (|82|i, i.e., Eq. (58 i with e = 0. However, 
they are the dynamics of a rare event when e ^ 0; because they are the most probable 
trajectory conditioned upon the rare event being observed. In fact, the H is the uni- 
form amount of excess kinetic energy required, added or reduced, to lead a rare event 
to occur. The Hamitonian dynamics, in this sense, has no forward predictive power 
based on given position and momentum; but it can predict detailed dynamical paths 
retrospectively, as in the Lagrangian formulation of classical mechanics and Fermat's 
principle for optics. 

It is easy to see from Eq. (83i that one conserved quantity in dynamics is q 2 — 
||6(<7)|| 2 . Hence, such a solution has a "uniform excess kinetic energy" compared with 
Eq. ( [82) 1. Conditioned on a given rare event, the problem of most probable ensem- 
ble is precisely the subject of Boltzmann-Gibbs' statistical mechanics and the Gibbs 
conditioning in the theory of large deviations ll39l[3Tl . 

Evolution of the landscape. The large-deviation rate function u st (x), as the sta- 
tionary solution to the HJE (|67|) or the "Boltzmann factor" for the stationary solution 



to Eq. ((59): 

u st (x) = lim--ln/ st (a;;e), 

e->0 e 

has been known for a long time to have a Lyapunov property for the dynamics dx/dt = 

j t u st (x(t)) = Vu st (x) ■ b(x) = -\\Vu st {x)f < 0. 

This means the ergodic, stationary stochastic dynamics contains a great deal of infor- 
mation on the time-dependent behavior of the system. In physics and biology, there 
is a growing interest to use the function u st (x) as an analytical visualization tool for 
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global behavior of complex dynamics: a landscape ROl [191 RTl l35l [T8l . u st (x) can 
even exist for vector fields b(x) which are non-conservative. The time-dependent HJE 
( |67| , therefore, can be interpreted as the evolution of the landscape. 
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Figure 1: The local minimum and its dynamics of u(9, t) on a circle, x is the location 
of a minimum and y its curvature, with their dynamics given in Eq. < (37] >. When the 
dynamical system on the circle has a stable fixed point at x\ i.e., b(p) = and 
b'(x*) < 0, the local minimum will approach x = x' and y = —b"(x'). The black 
filled circle in (A) represent a stable fixed point of the xy planar dynamics, the grey 
and open circles represent saddles and unstable node. However, when there is Hopf 
bifurcation as shown in (B), the x(t) and y(t) will be periodic functions of t. u(9, t) 
will not reach stationarity. However, as shown in Eq. (40 1, infinite time u(0,i) is a 
trival periodic function with a constant value. 
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Figure 2: Phase portrait for Hamiltonian system (42 1: q = f(p,q),p = g(p,q). The 
nullcline for g(p, q) = are p = and q = q* where q*s are where db(q)/dq = 0, 
shown in blue dashed lines. The nullcline for f(p, q) = are p — —b(q)/2, shown 
in red dashed lines. The orange lines of p — and p = —b(q) are the contour for 
H(q,p) = p 2 + b(q)p = 0. The intersections between blue and red dashed lines are 
equilibrium points of the Hamiltonian system: Centers are by green open circles and 
saddles are in black filled circles. Trajectories in green have H < and in pink have 
H > 0. 
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